Method of and an apparatus for processing seismic data

ABSTRACT

A method of and an apparatus for processing seismic data A method is provided of processing multi-component seismic data acquired by emitting multi-component seismic energy at a source location and acquiring seismic data at a multi-component seismic receiver located at a greater depth than the source location. The seismic data are decomposed into up-going constituents and down-going constituents ( 52 ). A de-signature and de-multiple operator is calculated ( 53 ) from the down-going constituents of the seismic data and from properties of the medium surrounding the receiver.

FIELD OF THE INVENTION

The present invention relates to a method of processing seismic data, inparticular multicomponent seismic data, to remove unwanted events fromthe acquired data. It also relates to an apparatus for processingseismic data. It also relates to a method of and an apparatus forcalculating a de-signature and de-multiple operator or a de-signatureoperator.

THE RELATED ART

FIG. 1(a) is a schematic illustration of the principles of a seismicsurvey. This seismic survey is intended to provide information about atarget geological reflector 3 disposed within the earth's interior.

The seismic survey shown in FIG. 1 includes a seismic source 4 disposedon the earth's surface 1. A seismic sensor 5, referred to hereinafter asa “receiver”, is also disposed on the earth's surface, at a distancefrom the seismic source 4. In use, the seismic source 4 is caused toactuate a pulse of seismic energy, and emitted seismic energy isdetected by the receiver 5.

FIG. 1(a) illustrates a land-based seismic survey. However, seismicsurveying of this general type is not restricted to land, and may becarried out in a marine environment or in the land-sea transition zone.For example, marine seismic surveying arrangements are known in whichone or more seismic sources are towed by a survey vessel; in such anarrangement, the receivers may be disposed on the sea-bed (a so-called“ocean bottom cable” survey or OBC survey), or the receivers may also betowed by a survey vessel. Furthermore, a land-based seismic survey isnot limited to the configuration shown in FIG. 1 (a), and it is possiblefor the seismic source or the seismic receiver to be disposed within theearth's surface. For example, in a vertical seismic profile (VSP)seismic survey a seismic source is placed on the earth's surface and areceiver is disposed within the earth's interior in a bore hole. In areverse VSP survey a seismic source is disposed within a bore hole, anda receiver is disposed on the earth's surface.

Only one seismic source 4 and one seismic receiver 5 are shown in FIG.1(a) for ease of explanation but, in general, a practical seismic surveywill contain an array of sources and an array of receivers.

One problem that arises in seismic surveying is that seismic energy maytravel from the source to the receiver along many paths. One reason forthis is that many other reflectors exist within the earth in addition tothe target reflector. In FIG. 1(a) this is illustrated schematically bya reflector 2 that overlies the target reflector 3. These additionalreflectors generate paths of seismic energy from the source to thereceiver that involve a reflection at a reflector other than the targetreflector. Another reason for the existence of multiple paths of seismicenergy is that seismic energy that is propagating upwardly within theearth will undergo reflection at the earth's surface 1, owing to thedifferent seismic properties of the earth and the air. This leads to theexistence of paths of seismic energy from the source to the receiverthat involve more than a single reflection at the target reflector.These paths give rise to unwanted events in seismic data acquired at thereceiver. An event in acquired seismic data that relates to seismicenergy that has undergone multiple reflections will hereinafter bereferred to a “multiple event”.

FIG. 1(a) illustrates the primary path of seismic energy for thissurveying arrangement, in which the path of seismic energy from thesource 4 to the receiver 5 involves only a single reflection, at thetarget reflector 3. (Refraction at the overlying reflector 2 has beenomitted from FIG. 1(a) for clarity.) Although seismic energy propagatingalong the primary path passes through the overlying reflector 2 on itsdownwards path from the source 4 to the target reflector 3, and again onits upwards path from the target reflector 3 to the receiver 5, theprimary path does not involve reflection at the overlying reflector 2.In an ideal seismic survey, only seismic energy that travelled along theprimary path would be detected by the receiver 5.

In a practical seismic survey, the receiver 5 will detect seismic energythat has travelled from the source 4 along many paths other than theprimary path. Examples of these other paths of seismic energy are shownin FIG. 1(b) to 1(d). FIGS. 1(b) to 1(d) illustrate paths of seismicenergy that involve more than one reflection, and these are known as“multiple events”. In FIG. 1l(b) downwardly propagating seismic energyfrom the source 4 is reflected by the overlying reflector 2 so that ittravels upwards to the earth's surface 1. The seismic energy is furtherreflected downwards at the earth's surface, and is then incident on thetarget reflector. An event of this general type is known as a“source-leg multiple”, since the additional reflections occur in thepath of seismic energy from the source to the target reflector.

FIG. 1(c) illustrates a seismic energy path in which the seismic energytravels direct from the source to the target reflector 3, and isreflected upwards at the target reflector 3. However, the reflectedseismic energy is not incident direct on the receiver, but is reflecteddownwards at the earth's surface, and is then reflected upwards at theoverlying reflector 2 before reaching the receiver. A seismic path ofthis type is known as a “receiver-leg multiple”, since the additionalreflections occur on the path of seismic energy from the targetreflector to the receiver.

FIG. 1(d) illustrates a seismic path in which seismic energy from thesource is incident on the target reflector 3, is reflected upwards tothe earth's surface, is reflected downwardly and undergoes a furtherreflection at the target reflector 3 before reaching the receiver. Inthis seismic energy path, the additional paths occur between the path ofseismic energy from the source to the target reflector and the path ofenergy from the target reflector to the receiver.

The seismic energy acquired at the receiver in a practical seismicsurvey will include events corresponding to the desired primary path1(a), but will also contain events relating to unwanted multiple pathssuch as the paths shown in FIGS. 1(b) to 1(d). In order to provideaccurate information about the target reflector it is desirable to beable to identify and remove multiple events from the seismic energyacquired at the receiver.

The seismic sources used in a land-based survey are normally vibrated orexplosive sources. If vibrators are used it is possible to perform amulti-component survey, using a multi-component vibrator that producesthree orthogonal source motions (two in orthogonal horizontal directionsand one in the vertical direction). If a seismic receiver is used thatcan record particle motion in three orthogonal directions, then it ispossible to perform a 3C×3C (or 9C) seismic survey. A suitable receiverfor this is one that can measure three orthogonal components of theparticle motion at the receiver, for example a receiver that containsthree orthogonal geophones—two geophones for measuring two orthogonalhorizontal components of the particle motion at the receiver, and athird geophone for measuring the vertical component of the particlemotion at the receiver.

A further problem in analysing the results of a seismic survey is thatmulti-component vibrators, which operate by imposing tractions on theearth's surface, emit three distinct wave types, known as P-, Sv- andSh-waves (P-waves are pressure waves, and Sv and Sh waves are shearwaves). The relative amplitudes of these different wave-types in seismicenergy emitted by a multi-component vibrator vary depending on thedirection of propagation of the seismic energy. A multi-componentreceiver records the three wave types, with a sensitivity that dependson the angle of incidence of the received seismic energy. When ageophone measures a component of the wavefield at the earth's surface,both P-waves and the two-types of S-waves are recorded withoutdistinction.

This is illustrated schematically in FIG. 4(a). FIG. 4(a) shows aland-based seismic survey in which a three-component seismic source 4—inthis case a multi-component vibrator—and a three-component seismicreceiver 5 are disposed on the earth's surface. As indicated in FIG.4(a), the seismic source 4 emits both P waves and S waves (only oneS-wave type is shown for clarity), and the receiver 5 detects bothP-waves and S-waves. The four views in FIG. 4(a) show, from left toright:

-   -   (i) Seismic energy generated by a horizontal source motion at        the vibrator 4, and being received by a horizontally-oriented        geophone at the receiver 5;    -   (ii) Seismic energy being generated by a horizontal source        motion at the vibrator 4, and being detected by a        vertically-oriented geophone at the receiver 5;    -   (iii) Seismic energy being generated by a vertical source motion        at the vibrator 4 and being detected by a horizontally-oriented        geophone at the receiver 5; and    -   (iv) Seismic energy being generated by a vertical source motion        at the vibrator 4 and being detected by a vertically-oriented        geophone at the receiver 5.

FIG. 4(a) illustrates only one horizontal component and one verticalcomponent of the source motion at the vibrator 4 and one horizontalgeophone component and one vertical geophone component at the receiver5. As noted above, in a full multi-component survey the vibrator 4 willfurther generate seismic energy by a source motion out of the plane ofthe paper, and the receiver will further comprise a third geophone thatdetects the particle motion along a line out of the plane of the paper.In total, thus, there are 9 combinations produced by the threeorthogonal source motions generated at the source and the threeorthogonal geophones at the receiver.

As noted above, many seismic receivers record P- and S-waves withoutdistinction, so that a seismic trace acquired at the receiver willincludes events due to received P-waves and events due to receivedS-waves. In many cases it is desirable to separate the P-events in aseismic trace from the S-events, since this provides additionalinformation about the earth's interior. In many cases, a geologicalstructure will have a different effect on P-waves than on S-waves. Theprocess of separating the P-events in a seismic trace from the S-eventsis generally referred to as decomposing the seismic trace into its P-and S-components.

There are a number of prior art approaches to decomposing seismic energyacquired at a receiver into P-components and S-components. There arealso a number of prior approaches to eliminating the effect of multiplereflections from the acquired seismic data.

For land-seismic data, C.P.A. Wapenaar et al, in “Decomposition ofmulticomponent seismic data into primary P- and S-wave responses”,Geophys. Prosp. Vol. 38, pp633-661 (1990) and P. Herrmann, in“Decomposition of multicomponent measurements into P- and S-waves”,Ph.D. thesis, Delft University of Technology (1992) have given anelastic decomposition scheme which decomposes data both at the receiverside (common shot gather) and at the source side (common receivergather). The wave field decomposition applied to a common source gatheralong the receiver line replaces the original particle velocitydetectors by pure P- and S- wave detectors. The wave field decompositionapplied to a common receiver gather similarly replaces the originalvibrator sources by pure P- and S-wave sources. Thus, the totaldecomposition scheme provides 9 different data sets: P, Sv and Shprocessed gathers from simulated P, Sv and Sh sources.

After the wave field decomposition, Wapenaar et al (1990, supra) applyan inverse scheme that eliminates the response of the earth's surfacefrom the decomposed results. In this technique, the earth's surface isassumed to be an interface between a solid and a vacuum. Eliminating theresponse of the earth's surface eliminates all multiple paths thatinvolve reflection at the earth's surface from the decomposed data, sogiving the primary P- and S-wave responses of the earth's interior.

This prior technique requires that the seismic sources are point sourcesand have a known wavelet or have a wavelet that can be estimated fromthe acquired data. The above prior art technique has been generalised byE. Holvik and L. Amundsen, in “decomposition of multi component seafloor data into primary PP, PS, SP, and SS wave responses”, ExpandedAbstracts of 68^(th) Annual Int Mtg of Society of ExplorationGeophysicists, pp2040-2043 (1998), to the case of a marine seismicsurvey in which ideal vibrators (traction sources) and geophones aredeployed at the sea floor. This technique relates to decomposing theacquired data, and removing events arising to reflections within thewater layer.

A further prior art technique has been developed by K. Matson and A.Weglein, in “Removing of elastic interface multiples from land and oceanbottom seismic data using inverse scattering” in Expanded Abstracts of66th Annual Int Mtg of Society of Acceleration of Geophysicists,pp1526-1529 (1996), and by K. Matson in “An inverse scattering seriesmethod for attenuating elastic multiples from multicomponent land andocean bottom seismic data”, Ph.D. thesis, University of British Columbia(1997). In this technique, inverse scattering series are used to developelastic schemes that attenuate multiple reflections from multi-componentland or ocean bottom seismic data.

The techniques developed by Holvik and Amundsen and by Matson andWeglein again require that the seismic sources must be point sourceswith a wavelet that is known or that can be estimated from the acquireddata.

L. Amundsen has proposed, in “Geophysics” Vol 66 pp 327-341 (2001), amethod of eliminating free-surface multiples from marine seismic dataacquired using a multicomponent seismic receiver disposed in the watercolumn or on the sea-bed. However, this method does not remove theeffects of multiple reflections associated with the sea-bed.

SUMMARY OF THE INVENTION

A first aspect of the present invention provides a method of processingmulti-component seismic data, the data having been acquired by emittingmulti-component seismic energy at a source location; and acquiringseismic data at a multi-component seismic receiver located at a greaterdepth than the source location, the method comprising the steps of:decomposing the acquired seismic data into up-going constituents anddown-going constituents; and calculating a de-signature and de-multipleoperator from the down-going constituents of the acquired seismic dataand from properties of the medium surrounding the receiver.

The de-signature and de-multiple operator of the invention is effectiveat attenuating or completely removing the effects of the overburden fromthe seismic data. It is effective at removing from the seismic data allmultiple reflections associated with any interface above the level ofthe receiver or with any interface. at the receiver level. (The term“interface” covers any discontinuity in acoustic or elastic propertiesthat causes partial reflection of seismic energy.)

It is also effective at attenuating or completely removing the effectsof the source radiation characteristics (or the source “signature”) fromthe data.

The method may comprise processing the acquired seismic data using thede-signature and de-multiple operator thereby to attenuate or removeseismic events arising from multiple reflections. This may be done byprocessing the up-going constituents of the acquired seismic data usingthe de-signature and de-multiple operator. Alternatively, thede-signature and de-multiple operator may be applied to the complete(that is, undecomposed) acquired seismic data.

The step of processing the acquired seismic data may comprise selectinga desired seismic signature for the source.

The method may further comprise decomposing the seismic data into P-waveand/or S-wave data. This decomposition step may be performed on thede-signatured, de-multipled data. Alternatively, the acquired data maybe decomposed into P-wave and/or S-wave data before the de-signature andde-multiple operator is applied.

The decomposition into P-wave and/or S-wave data may be a receiver-sidedecomposition and/or a source-side decomposition.

A second aspect of the invention provides a method of processingmulti-component seismic data, the data having been acquired by emittingmulti-component seismic energy at a source location; and acquiringseismic data at a multi-component seismic receiver located at a greaterdepth than the source location, the method comprising the steps of:

-   -   decomposing the acquired seismic data into up-going constituents        and down-going constituents; and    -   calculating a de-signature operator from the initial down-going        constituents of the acquired seismic data and from properties of        the medium surrounding the receiver.

The initial down-going constituents of the seismic data acquired at thereceiver relate to seismic energy that has travelled direct from thesource to the receiver without undergoing any reflection. They maytherefore be used to calculate a de-signature operator for attenuatingor eliminating the effects of the radiation characteristics of theseismic source on the acquired data.

The methods of the invention may be applied to pre-existing seismicdata. Alternatively, the method may further comprises the steps of:emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location.

Other aspects of the invention are claimed in independent claims 10, 12,13, 14, 16, 23, 24, 26, 27,28 and 30.

BRIEF DESCRIPTION OF THE DRAWINGS

Preferred embodiments of the present invention will now be described byway of illustrative example with reference to the accompanying figuresin which:

FIG. 1(a) is a schematic view of a seismic survey showing the desiredpath of seismic energy;

FIGS. 1(b), 1(c) and 1(d) are schematic views of a seismic surveyshowing three undesired paths of seismic energy;

FIG. 2 is a flow diagram showing possible steps of a method of obtaininga de-signature and de-multiple operator according to an embodiment ofthe present invention;

FIG. 3(a) is a schematic view of a real seismic survey;

FIG. 3(b) shows a hypothetical seismic survey corresponding to theseismic survey of FIG. 3(a);

FIG. 3(c) shows a further hypothetical seismic survey that is areciprocal of the hypothetical seismic survey shown in FIG. 3(b);

FIGS. 4(a) to 4(d) illustrate the effect of applying a method of thepresent invention to a land seismic survey;

FIGS. 5(a) and 5(b) illustrates the effects of applying the method ofthe present invention to an ocean bottom seismic survey;

FIGS. 6(a) and 6(b) illustrate the effect of applying a method of thepresent invention to a land seismic survey in which the receiver isdeployed within the earth's interior;

FIGS. 7(a) to 7(d) show raw synthetic seismic data;

FIGS. 8(a) to 8(d) show the up-going component of the synthetic seismicdata of FIGS. 7(a) to 7(d);

FIGS. 9(a) to 9(d) show the down-going part of the synthetic seismicdata of FIGS. 7(a) to 7(d);

FIG. 10(a) to 10(d) show the particle velocities obtained from thesynthetic seismic data of FIGS. 7(a) to 7(d) using a method of thepresent invention;

FIG. 11(a) to 11(d) show the PP, SP, PS and SS wavefields obtained fromthe synthetic seismic data of FIG. 7(a) to 7(d) using a method of thepresent invention;

FIG. 12 is block schematic view of an apparatus according to theinvention.

BRIEF DESCRIPTION OF THE PREFERRED EMBODIMENTS

The principle of one method of obtaining a de-signature and de-multipleoperator of the present invention is illustrated in FIGS. 3(a) to 3(c).FIG. 3(a) illustrates an actual seismic survey. This seismic surveycomprises a receiver x_(r) (14) that is located at a depth z_(r). Themedium within which the seismic survey is carried out is classified intothe overburden for all depths less than the receiver depth (that is forall z<z_(r)), and as the sub-surface for depths greater than the depthof the receiver (that is, for all z>z_(r)). The overburden is anarbitrary inhomogenous medium.

A seismic source x_(s) (13) is disposed in the overburden, at a depthless than the depth of the receiver. The seismic source emits bothup-going and down-going seismic energy.

The inhomogenous nature of the overburden is indicated schematically bythe presence of interfaces 10, 11, 12 within the overburden. The seismicproperties of the overburden change at these interfaces. In the case ofa marine seismic survey, the upper interface 12 might be the sea-surface(that is, an air/water interface). Because of the existence of theseinterfaces in the overburden, up-going seismic energy is reflecteddownwards in the overburden so that multiple events occur in the seismicdata acquired at the receiver 14.

The seismic source x_(s) (13) emits both S- and P-waves. The receiverx_(r) (14) records receive S- and P-waves without distinction. Thereceiver may be any device that measures particle motion, such asparticle displacement, particle velocity (the time derivative ofparticle displacement), or particle acceleration (the time derivative ofparticle acceleration). One commonly-used seismic receiver is ageophone, but the invention is not limited to processing seismic dataacquired using one or more geophones as the receiver(s).

FIG. 3(b) illustrates an ideal seismic surveying arrangement thatcorresponds to FIG. 3(a). In the idealised arrangement of FIG. 3(b), theoverburden is an homogenous medium, so that up-going seismic energy fromthe seismic source continues to propagate upwards and is not reflecteddownwards, so that seismic energy acquired at the receiver does notcontain multiple events. The interfaces 10, 11 and 12 in the overburdenof FIG. 3(a) have been illustrated in broken lines in FIG. 3(b) toemphasise that the overburden is homogenous in FIG. 3(b).

In the ideal seismic surveying arrangement of FIG. 3(b) the seismicsource 13 is centred at a location x_(r) ⁻, which is at a depth z_(r) ⁻.The depth of the source, z_(r) ⁻, is at an infinitesimal distance Eabove the depth z_(r) ⁻ at which the receiver 14 is positioned. Thedistance E is vanishingly small.

FIG. 3(c) shows an alternative idealised seismic surveying arrangement.This corresponds to the idealised seismic surveying arrangement of FIG.3(b), except that the positions of the seismic source 13 and the seismicreceiver 14 have been interchanged. That is, in FIG. 3(c) the seismicsource 13 is located at position x_(r), at depth z_(r), and the seismicreceiver 14 is located at location x_(r) ⁻ at a depth z_(r) ⁻ which is adistance E above the depth z_(r) of the source.

Since the seismic surveying arrangement of FIG. 3(c) is identical to theseismic surveying arrangement of FIG. 3(b) except that the positions ofthe source and receiver have been interchanged, the theorem ofreciprocity requires that the seismic surveying arrangement shown inFIG. 3(b) will produce the same seismic data as the seismic surveyingarrangement shown in FIG. 3(c).

In one embodiment of the invention, Betti's reciprocity theorem is usedto transform actual seismic data required in a real seismic survey ofthe type shown in FIG. 3(a) into the seismic data that would be obtainedin an idealised seismic survey of the type shown in FIG. 3(b). Thisenables the effect of the inhomogenous overburden to be eliminated,thereby eliminating multiple events.

The present invention also provides decomposition operators, thatdecompose the acquired seismic data into PP, PS, SP, and SS components.

Betti's reciprocity theorem gives an integral equation relationshipbetween two independent elastic wavefields that are defined in aspecific volume enclosed by a mathematical or physical surface. Therelationship between the two wavefields is governed by possibledifferences in the parameters of the medium, by possible differences inthe source distribution, and by possible differences in boundaryconditions. In the invention, Betti's theorem is used to provide anintegral equation procedure that can be used to transform a wavefieldthat was recorded in an actual seismic survey with an overburdenresponse (that is, a survey in which multiple reflections occur owing toan inhomogenous overburden) into wavefields that would have beenrecorded in an idealised hypothetical seismic experiment with theoverburden response absence (the overburden in the hypothetical seismicsurveying arrangement is homogenous and so produces no response). Thepresent invention makes it possible to remove undesired events in theseismic data that arise owing to multiple reflections of the types shownin FIGS. 1(b) to 1(d), so that a multiple-free image can be obtained.The present invention enables all multiple reflections associated withany interface above or at the level of the receiver to be removed fromthe seismic data.

Other than the position of the source, or sources, the present inventiondoes not require any source characteristics for the calculation toeliminate the waves scattered from the overburden. The radiationcharacteristics of the physical seismic source(s) is eliminated in thetransformation from the physical seismic surveying arrangement to theidealised hypothetical seismic surveying arrangement. The source in thehypothetical idealised surveying arrangement is assumed to be a singlepoint source that is generates components orthogonally oriented withrespect to one another, and has a desired wavelet.

This embodiment of the method of eliminating the response of theoverburden and the effect of the multi-component source radiationpattern from the acquired seismic data is referred to hereinafter as“Betti designature and elastic demultiple”. A Betti designature andelastic demultiple method of the invention has the following advantages:

-   -   (a) It preserves the amplitude of the primary reflection while        eliminating all events in the seismic data that arise owing to        scattering of up-going energy in the overburden;    -   (b) It requires no knowledge of the medium below the depth of        the receiver;    -   (c) It requires no knowledge of the medium above the depth of        the receiver; and    -   (d) It requires information only about the local density and        local velocities of elastic wave propagation at the depth of the        receiver.

In a method of the present invention recorded seismic data aretransformed to new, idealised data that would be recorded in ahypothetical seismic survey in which the overburden response was absent(that is in which the overburden was completely homogenous). The sourcein this hypothetical seismic survey is a point source, with some desiredsignature radiation characteristics (wavelet). In the processed data,the effect of the physical source and its radiation characteristics isabsent. In other words, the new data have been designatured.

A method of the present invention may be applied to data obtained in anyseismic survey that uses a multi-component seismic source and amulti-component seismic receiver. Three particular seismic surveyingarrangements to which the invention can be applied are as follows:

A. A 3C×3C (or 9C) land seismic survey, in which the seismic sourcegenerates three orthogonal source motions (usually two horizontal sourcemotions and one vertical source motion). The seismic wavefield isrecorded by a receiver that contains three orthogonal geophones,deployed at or below the earth's surface. The three orthogonal geophonesmeasure three orthogonal components of the particle velocity vector(usually the x-, y- and z-components).

B. A 3C×4C (or 12C) marine seismic survey, in which the seismic sourcegenerates three orthogonal source motions (usually two horizontal sourcemotions and one vertical source motion). The seismic wavefield isrecorded by a receiver that contains three orthogonal geophones,deployed just below the seabed and one hydrophone deployed just abovethe seabed. The three orthogonal geophones measure three orthogonalcomponents of the particle velocity vector, and the hydrophone measuresthe pressure field (which is a scalar quantity).

C. A 3C×6C (18C) bore hole seismic experiment in which the seismicsource again generates three separate orthogonal source motions. Thesource may be disposed on land, on the sea floor, or in the watercolumn. The seismic receiver preferably measures three orthogonalcomponents of the particle velocity, and three orthogonal components ofthe vertical traction vector or estimates thereof.

According to the invention, the de-multiple and de-signature operator isobtained from the down-going constituents of the acquired seismic data.This means that the measurements made at the receiver must be sufficientto allow the acquired wavefield to be decomposed at the receiver sideinto up-going and down-going wave constituents at the receiver. In themost general case, a proper wave field decomposition requires that thereceiver measures the three components of the particle velocity vectoras well as the three components of the vertical traction vector. Forland seismic data, however, the vertical traction vector is zero, sothat only the three components of the particle velocity vector need tobe measured.

For ocean bottom seismic data, the two horizontal components of thevertical traction vector are zero, and the vertical component of thevertical traction vector just below the sea-bed is equal in magnitudeand opposite in sign to the pressure field just above the sea floor.Thus, providing a hydrophone just above the sea-bed enables the verticalcomponent of the vertical traction vector (which is the only non-zerocomponent of the traction vector) to be determined.

Applying the present invention to a bore hole seismic survey requires,in principle, that the three components of the vertical traction vectormust be measured or otherwise be known. In most cases, however, thethree components of the vertical traction vector are not known. In orderto overcome this, it is usual to record the three-component particlevelocity at two different depths within the bore hole. Thus, thereceiver system is essentially a six-component receiver system, so thata bore hole seismic survey may be considered to be an 18C survey.

FIG. 2 is a flow diagram showing the principal steps of a methodaccording to the present invention.

The invention may be carried out on suitable pre-existing seismic data.In this case, the method would commence with the step (not shown) ofretrieving suitable seismic data from storage. Alternatively, the methodmay commence with the step (not shown) of performing a multi-componentsource, multi-component receiver seismic survey to acquire suitableseismic data.

Initially at step S1 the local elastic properties of the medium in thevicinity of the receiver are determined. For a marine seismic survey, asone example, this step requires a determination of the local elasticproperties of the sea-bed in the vicinity of the receiver. Theseproperties may be determined from measurements, or they may be estimatedfrom knowledge of the geological structure of the survey site.

At step S2, the multi-component seismic data that has been acquired orretrieved from storage is decomposed into up-going and down-goingwavefield constituents. In principle step S2 could be carried outbefore, or simultaneously with, step S1. It is generally preferable toperform step Si before step S2, however, since the local elasticproperties obtained from step Si may then be used to decompose theseismic data into up-going and down-going waves.

At step S3 the de-signature and de-multiple operator is calculated fromthe multi-component data and from the elastic properties of the medium.

At step S4 a desired signature for the seismic source is selected. Thisdesired source signature, the multi-component data obtained at S2, andthe operator calculated at S3 are then used to solve an integralequation to find the de-signatured and de-multipled wavefield, at stepS5. In essence, S5 consists of transforming the seismic data obtained ina real seismic survey of the type shown in FIG. 3(a) into the data thatwould be obtained in an idealised seismic survey shown in FIG. 3(b). Ina preferred embodiment Befti's theorem is used to transform the actual,measured wavefield enclosed by a surface S of the seismic surveyingarrangement of FIG. 3(a) into the wavefield enclosed by the same surfaceS in the idealised seismic surveying arrangement of FIG. 3(b). Thesurface S is made up of the surfaces Σ and S_(R) shown in FIGS. 3(a) to3(c). The surface Σ is a horizontally plane surface that is located atdepth z_(r) ⁻—that is, located an infinitesimal distance above the depthof the receiver (which is located at a depth z_(r)). The surface S^(R)is a hemispherical surface having a radius of R.

In Step S5 the de-signature and de-multiple operator may be applied tothe up-going wavefield constituents obtained at step S1. Alternatively,the de-signature and de-multiple operator may be applied to the seismicdata as originally-acquired at the receiver.

Finally, at step S6, the de-signatured and de-multipled data isdecomposed into the primary PP, PS, SP, and SS wave responses that wouldbe recorded if pure P- and S-wave sources and pure P- and S-wavereceivers had been used in the seismic survey.

In one embodiment, Step S6 has two principal operations. Initially, themulti-component data are Fourier transformed over horizontal receiverco-ordinates. The transformed data then undergo receiver side wavefields decomposition, in the wave number domain, into up-going pressureand shear wave fields. This is done by operating on the data by areceiver decomposition matrix. The decomposed data then undergoes areverse Fourier transform, to give the up-going pressure and shear wavefields measured at the receiver.

The thus-obtained multi-component data are then Fourier transformed overhorizontal source co-ordinates. The transformed data then undergoessource side wave field decomposition in the wave number domain, intodown-going pressure and shear wave fields. This is done by operating onthe data with a source decomposition matrix. The data then undergoes areverse Fourier transform, to give the down-going pressure and shearwave fields emitted at the source.

The combined elastic wavefield decomposition on the source and receiversides gives seismic data that are equivalent to data that would beobtained in an idealised survey in which the overburden response wasabsent, and in which single component P- or S-wave sources and singlecomponent P- or S-wave receivers were used.

A further advantage of the invention is that it removes, or at leastsignificantly attenuates, ground-roll noise (in the case of a land-basedseismic survey) or Scholte waves (in the case of a marine seismicsurvey). Ground roll noise is often the largest source of noise in aland-based seismic survey, so the removal of the effects of ground rollnoise is highly desirable.

Step S6 can be omitted if desired. Alternatively, the step S6 ofdecomposing the data into the PP, PS, SP, and SS wavefields may becarried out before the step S5 of applying the de-signature andde-multiple operator. If the decomposition into PP, PS, SP and SSwavefields is done first, this step would replace the step S2 ofdecomposing the seismic data into up-going and down-going waves is notrequired.

If the medium in which the seismic survey is carried out is horizontallylayered, the Betti designature and elastic demultiple scheme, followedby the elastic source-receiver decomposition scheme, is greatlysimplified. It may be conveniently implemented in the tau-p domain or inthe frequency-wave number domain.

The method is preferably carried out on common shot gathers, or oncommon receiver gathers when source array variations are negligible. Inthe latter domain the assumption of a layered earth may possibly degradethe results, but in most cases this is not very significant.

Once the seismic data has been corrected for the effects of theoverburden and the source radiation characteristics using thede-signature and de-multiple method of the invention, and optionally hasbeen decomposed into P-wave and/or S-wave constituents, it may then besubjected to further processing steps such as, for example, one or moreconventional processing steps.

FIG. 4 illustrates the effect of the present invention on a land seismicsurvey in which a multi component seismic source 4 and a multi-componentseismic receiver 5 are disposed on the earth's surface.

FIG. 4(a) illustrates a real seismic survey. When the source 4 isactuated, down-going P- and S-waves are emitted. P-waves are denoted bysolid lines, and S-waves are denoted by broken lines. The arrow in thebox denoting the source (labelled “src”) indicates the direction ofsource motion generated by the source. The seismic energy emitted fromthe source 4 undergoes reflection within the earth, and the up-goingreflected waves are incident on the receiver 5. The receiver 5 containsat least two orthogonal geophones, and the arrow in the box denoting thereceiver (labelled “rec”) indicates the component of particle motionthat is measured by the receiver. When the receiver records thewavefield at the earth's surface it records both p-waves and s-waveswithout distinction.

FIG. 4(a) shows four cases

-   -   (i) the source generates a horizontal particle motion, and the        receiver records particle motion in that horizontal direction;    -   (ii) the source generates a horizontal particle motion and the        receiver records the vertical component of the received particle        motion;    -   (iii) the source generates a vertical particle motion, and the        receiver records a horizontal component of the received particle        motion; and    -   (iv) the source generates a vertical particle motion, and the        receiver records the vertical component of the received particle        motion.

FIG. 4(b) illustrates the results of applying the de-signature andelastic de-multiple process to the data recorded in FIG. 4(a). That is,FIG. 4(b) shows the results after the seismic data acquired in theseismic surveying arrangement of FIG. 4(a) has been processed accordingto steps S1 to S5 of FIG. 2. The effect of this processing is toeliminate the response of the overburden (that is, the response of theearth's surface in this case) from the acquired seismic data. Theprocessing also eliminates the effect of the radiation characteristicsof the seismic source 5. The source is now a point source of seismicenergy. In principle energy may propagate upwards past the source andreceiver, since the overburden has been replaced by an infinitehomogenous medium, and this is indicated in FIG. 4(b).

As noted above, the source 4 emits both P- and S-waves. In FIGS. 4(a)and 4(b), the receiver 5 records both S-waves and P-waves withoutdistinction.

FIG. 4(c) illustrates the results of applying wave fields decompositionat the receiver side to the seismic data of the idealised hypotheticalseismic surveying arrangement shown in FIG. 4(b). This has the effect ofreplacing the geophones used as the receiver in FIGS. 4(a) and 4(b) by asimulated pure P-wave detector (FIGS. 4(c)(i) and 4(c)(iii)), or by asimulated pure S-wave detector (FIGS. 4(c)(ii) and 4(c)(iv)).

Finally, FIG. 4(d) shows the results of further applying a source sidewavefield decomposition. This has the effect of replacing the seismicsource used in FIGS. 4(a) to 4(d), which emits both P-waves and S-waves,by a simulated pure P-wave source

(FIGS. 4(c)(i) and 4(c)(ii))or by a simulated pure S-wave source (FIGS.4(c)(ii) and 4(c)(iv)).

FIG. 4(d) illustrates the results of applying steps S1 to S6 of theembodiment shown in FIG. 2. As can be seen, the method has eliminatedthe overburden response, eliminated the radiation signature of thesource, and has decomposed the acquired data into PP, PS, SP, and SSevents.

FIG. 5 illustrates the results of applying a method of the presentinvention to an ocean bottom seismic survey. In the surveyingarrangement shown in FIG. 5(a) a seismic source 4 and a seismic receiver5 are each disposed on the sea-bed. The arrow in the box denoting thesource indicates the direction of source motion generated by the source.The receiver 5 contains at least two orthogonal geophones, and the arrowin the box denoting the receiver indicates the component of particlemotion that is measured by the geophone. (A receiver used in a practicalocean bottom seismic survey may also contain a hydrophone, but this isnot shown in FIG. 5.) The four combinations of horizontal and verticaldirection of source motion and horizontal and vertical component ofparticle motion measured at the geophone shown in FIGS. 5(a)(i) to5(a)(iv) correspond to the four cases shown in FIGS. 4(a)(i) to4(a)(iv).

As can be seen in FIG. 5(a), the source emits down-going S-waves andP-waves, which undergo reflection within the earth so as to be incidentupon the receiver 5. FIG. 5(a) also illustrates one possible water layermultiple path for seismic energy, in which upgoing P-waves are reflectedat the water surface so as to be incident on the receiver 5.

The geophones in the receiver 5 measures the incident S-waves andP-waves without distinction.

FIG. 5(b) illustrates the effect of applying a de-signature and elasticde-multiple method of the present invention to the seismic data acquiredin the seismic surveying arrangement shown in FIG. 5(a). As can be seen,the water layer multiple paths are eliminated, since the inventioneliminates the response of the overburden (that is, the water layer inthis case) from the acquired setsmic data. The method of the inventionalso eliminates the effect of the radiation signature of the source, andthe source in the idealised seismic surveying arrangement of FIG. 5(a)is a point source of force.

FIG. 5(b) illustrates the effects of applying steps S1 to S5 shown inFIG. 1 to the seismic data acquired in the seismic surveying arrangementof FIG. 5(a). The effect of applying further steps of the invention—thatis, applying the receiver side decomposition and the source sidedecomposition would provide identical results to those shown in FIGS.4(c)(i) to 4(d)(iv).

FIG. 6(a) and 6(b) illustrate the results of applying the presentinvention to a land seismic survey which the receiver and source aredeployed within the earth. In the seismic surveying arrangement shown inFIG. 6(a) the seismic source 4 and the seismic receiver 5 are bothdeployed within the earth, with the receiver being at a greater depththan the source. FIG. 6(a) illustrates the primary path of seismicenergy from the source to the receiver, in which down-going S-waves andP-waves are reflected by a target reflector 3 so as to be incident onthe receiver 5. FIG. 6(a) also shows unwanted multiple reflection pathsin which up-going S-waves and P-waves undergo reflection in layersdisposed above the source so as to be incident on the receiver.

In FIG. 6(a) the arrows in the box denoting the source 4 indicate thedirection of source motion generated by the source. The arrows in thebox indicating the receiver indicate the component of particle motionthat is measured by a geophone contained in the receiver. As in FIGS. 4and 5, the four combinations of horizontal and vertical source motionand horizontal and vertical geophone orientation are shown in FIG.6(a)(i) to 6(a)(iv).

FIG. 6(b) illustrates the effect of applying the present invention tothe seismic data acquired in the seismic survey arrangement of FIG.6(a). As in FIGS. 4 and 5, the effect of applying a de-signature andde-multiple method of the invention to the data obtained in the seismicsurvey arrangement of FIG. 6(a) is to eliminate the effect of theoverburden from the acquired seismic data. The effect of the radiationsignature of the seismic source 5 is also elirninated. Thus, FIG. 6(b)illustrates an idealised seismic survey in which the overburden responseis eliminated (corresponding to an homogenous overburden), and in whichthe seismic source 4 is a point source of force. FIG. 6(b) illustratesthe results of applying steps S1 to S5 of FIG. 2 to the seismic dataacquired in the real seismic surveying arrangement of FIG. 6(a).

The effect of applying a receiver side decomposition and a sourcedecomposition would correspond to the results shown in FIGS. 4(c)(i) to4(d)(iv) respectively.

FIGS. 7(a) to 11(d) illustrate results of applying the present inventionto simulated seismic data.

FIGS. 7(a) to 7(d) shows raw simulated seismic data. This data wassimulated for a land seismic survey in which the source and receiver areboth located on the earth's surface. The seismic data was simulatedusing a reflectivity code of the type proposed by B. L. Kennett in“Seismic Wave Propagation in Stratefied Media” Cambridge UniversityPress (1993). The model used to simulate the seismic data consists of aflat free surface on which the source and receiver are located. Belowthe free surface is a 300 m thick layer, having a flat lower interface.This layer is an isotropic elastic layer, with a P-wave velocity of 1800m/s, and s-wave velocity of 800 m/s, and a density of 200 kg/m³. Belowthis layer is an isotropic elastic infinite half space, having a p-wavevelocity of 2500 m/s, an S-wave velocity of 1500 m/s, and a density of2200 kg/m³. Since the model is horizontally layered, it is necessary toconsider a vertical source motion and one horizontal source motion, andto consider a geophone that records either a vertical component or ahorizontal component of the particle motion. Thus, it is only necessaryto consider a 2C×2C (4C) case, rather than the fall 3C×3C (9C) case.

In FIGS. 7(a) to 10(d), the notation Vz and VX indicate verticalparticle motion and horizontal particle motion at the receiver,respectfully. The terms F, and F, indicate that the geophone records thevertical component and a horizontal component of the particle motion,respectively.

FIG. 8(a) to 8(d) illustrates the up-going constituent of the raw dataof FIGS. 7(a) to 7(d) respectively. FIGS. 9(a) to 9(d) show thedown-going constituent of the raw seismic data of FIGS. 7(a) to 7(d)respectfully. FIGS. 8(a) to 9(d) were obtained by applying equations(A-28), (A-29), (A-32) and (A-33) below. The effect of the free surfaceis still present in FIGS. 8(a) to 9(d).

FIGS. 10(a) to 10(d) illustrate the effect on the seismic data of FIGS.7(a) to 7(d) of removing the effects of the free surface and the sourcesignature. This is done by applying a de-signature and de-multiplemethod of the present invention. That is, FIGS. 10(a) to 10(d) shows theeffect of applying steps Sl to S5 of FIG. 2 to the raw seismic datashown in FIGS. 7(a) to 7(d), respectively. It can be seen that seismicevents in FIGS. 7(a) to 7(d) that relate to multiple reflections havebeen eliminated from the seismic data shown in FIGS. 1 0(a) to 10(d).

The events in the seismic data of FIGS. 1 0(a) to 1 0(d) include all PP,PS, SP and SS events, as well as related head-waves.

FIGS. 1(a) to 11(d) illustrate the effect of applying a receiver-sidedecomposition and a source-side decomposition to the data shown in FIGS.1 0(a) to 1 0(d), respectively. The effect of this is to decompose theseismic data into PP, SP, PS and SS events. The seismic data shown ineach of FIGS. 11(a) to 11(d) includes only one event, together with therelated head-waves.

The present invention has been described above with reference to acombined de-multiple & de-signature method. In an alternativeembodiment, the invention provides a de-signature method only. In thisalternative embodiment, the operator is computed from the down-goingconstituent of the initial arrival of seismic energy at the receiver.The initial arrival constitutes seismic energy that has travelled directfrom the source to the receiver without undergoing any reflection, forexample along the path 6 shown in FIG. 1(a). The incident part of thewavefield may be determined by, for example, the method proposed foracoustic or elastic earth by A. Weglein and B. G. Secrest in “Waveletestimation for a multi-dimensional acoustic or elastic earth”,Geophysics Vol. 55, pp902-913 (1990).

Furthermore, if the signature/wavelet of the multicomponent source isknown or can be estimated, the present invention may also be used toprovide a de-multiple method only. The principal steps of such a methodare to (a) obtain a combined de-signature and de-multiple operatoraccording to a method of the invention, (b) obtain a de-signatureoperator from the known or estimated radiation characteristics of theseismic source and (c) obtain a de-multiple operator that will combinewith the de-signature operator obtained at (b) to produce the combinedde-signature and de-multiple operator obtained at (a).

A detailed mathematical description of one embodiment of the presentinvention will now be given.

Consider a volume V enclosed by the surface S=Σ+S_(R) withoutward-pointing normal vector n, as depicted in FIGS. 3(a) to 3(c). Σis a horizontal plane surface located at depth z_(r) ⁻ infinitesimallyabove the multicomponent receiver x_(r) 14 located at depth level z_(r).To simplify the analysis, it is assumed that the solid medium ishomogeneous and isotropic at depth z_(r) and in a infinitesimally thinregion below depth z_(r). The invention does not require knowledge ofthe properties of the overburden (at depths z<z_(r)) and the subsurface(at depths z>z_(r)), and the overburden and sub-surface may both bearbitrarily inhomogeneous, anisotropic and anelastic media. When seismicdata is acquired close to a physical surface, the receiver(s) is/arealways infinitesimally below the physical surface. Physical surfacestypically are the surface of the earth (which is assumed to be asolid/vacuum surface), the ocean bottom (which is assumed to be asolid/fluid boundary), or another interface in the earth (assumed to bea solid/solid surface). S_(R) is a hemispherical surface of radius R.The Cartesian coordinates will be denoted by x=(ξ, x₃), where ξ=(x₁,x₂). For notational convenience, x₃=z will also be used. The x₃-axis,which is positive downwards, is parallel to n. The x₁, x₂-axes are inthe Σ plane.

Initially, an integral relationship is established between themulticomponent source, multicomponent receiver data acquired in a realseismic survey—these data contain the scattering response of the mediumabove the receivers—and the desired multicomponent source,multicomponent receiver data that would be obtained in an idealisedseismic survey with that scattering response absent.

The seismic source used in the survey is assumed to separately generatethree orthogonal source motions. The desired multicomponent data arethose data that would be recorded in an idealised, hypotheticalmulticomponent seismic survey using three orthogonal point forces actingseparately, with a desired, equal signature when the medium above thereceivers is homogeneous and extends upwards to infmity, and hasparameters equal to those of the medium at the receiver depth level. Theoverburden in the idealised survey thus is an elastic, isotropichalfspace. The geology below the receiver level is naturally the same inthe physical and hypothetical seismic surveys. TABLE 1 State P State HState Ĥ Force component f_(n)(x_(s)) ãδ(x − x_(r) ⁻)δ_(in) ãδ(x −x_(r))δ_(im) Traction source 0 0 0 Particle velocity v_(mn)(x_(r)|x_(s)){tilde over (v)}_(mn)(x_(r)|x_(r) ⁻) {tilde over (v)}_(mn)(x_(r)⁻|x_(r)) Vertical traction s_(mn)(x_(r)|x_(s)) {tilde over(s)}_(mn)(x_(r)|x_(r) ⁻) {tilde over (s)}_(mn)(x_(r) ⁻|x_(r))Notation for sources and fields entering the elastodynatnic equationsrelated to states P, H and Ĥ depicted in FIGS. 3(a) to 3(c),respectively.

In the physical survey, the seismic source is unidirectional, either avibrator or a force. A vibrator source consists of one vertical vibrator(source of tensile stress) and two horizontally oriented vibrators(sources of shearing stress). The vibrator, located at position x₅ needsnot be an ideal point source. Its source signature need not be known. Inthe hypothetical surveys H and Ĥ, the sources are point sources of forcein the n and m direction, respectively, with source signature ã.

Consider the seismic survey with configuration as depicted in FIG. 3(a),which will be referred to as a “physical” seismic survey forconvenience. The recorded m:th component of particle velocity vector atreceiver location x_(r), infinitesimally below Σ, due to a source actingin direction n at centre coordinate x_(s) with unknown source strengthand radiation pattern is denoted by v_(mn). Likewise, the m:th componentof vertical traction vector is denoted by s_(mn). The source and fieldvariables for the physical seismic survey, which is denoted for short“state P”, are listed in Table 1. Note that the surface Σ may or may notcoincide with a physical surface. In the case that recording takes placejust below the surface of the earth, or just below the ocean bottom, thesurface Σ obviously coincides with a physical surface.

The desired wavefields, {tilde over (v)}_(mn) and {tilde over (s)}_(mn)that it is proposed to solve for are the responses of the medium fromthree orthogonal point forces with desired signature ã when the mediumabove the receiver level is a halfspace as shown in FIG. 3(b). In thiscase, the surface Σ is a non-physical boundary. The desired particlevelocity and vertical traction vector responses are recorded at locationx_(r) infinitesimally below Σ for the point forces located at x_(r) ⁻ onΣ. The source and field variables for this hypothetical seismicexperiment which is denoted, for short, “state H”, are listed in Table1.

To establish the integral relationship between the physical state P andhypothetical state H seismic experiments the hypothetical “state Ĥ” isintroduced. This is a further idealised seismic survey, and correspondsto the idealised seismic survey H if FIG. 3(b) except that the source 5and the receiver 6 are interchanged with one another. Thus, thewavefields {circumflex over (v)}_(mn) and ŝ_(mn), in the hypothetical“state Ĥ” are the reciprocal wavefields to the ones in state H, and obeythe reciprocity relation{circumflex over (v)} _(nm)(x _(r) ⁻ |x _(r))={tilde over (v)} _(mn)(x_(r) |x _(r) ⁻),   (1)ŝ _(nm)(x _(r) ⁻ |x _(r))={tilde over (s)} _(mn)(x _(r) |x _(r) ⁻),  (2)

Thus, {circumflex over (v)}_(mn) and ŝ_(mn) are responses at locationx_(r) ⁻ on surface Σ due to a point force, with signature ã, acting indirection m at location x_(r) infinitesimally below Σ as depicted inFIG. 3(c). The surface Σ is, as in the desired idealised state H, anartificial, non-physical boundary. Below, we shall use the fact thatseismic data acquired in a survey in the hypothetical state Ĥ consist ofupgoing events only scattered from the subsurface below Σ. In addition,the direct wavemodes from the sources to the receivers are upgoingevents since the sources are ifinitesimally below the receivers. Thesource and field variables for state Ĥ are listed in Table 1.

Next Betti's reciprocity theorem is applied in the volume V enclosed bythe surface Σ+S_(R→∞), where S_(R→∞) is located at infinity. Similarlyto, for instance L. Amundsen et al. in “Elimination of free-surfacerelated multiples without the need of the source wavelet”, Geophysics,Vol 66, pages 327-341 (2000), a frequency-domain integral equationdescribing the relationship between state P and state H is obtained:ãv _(mn)(x _(r) |x _(s))=∫dS(ξ)└ŝ _(im)(x|x _(r))v _(in)(x|x_(s))−{circumflex over (v)}_(im)(x|x _(r))s _(in)(x|x _(s))┘.   (3)

The surface S_(R) does not contribute to the integral as R→∞ (this isthe radiation condition, see Y. H. Pao and V. Varatharajulu in “Huygen'sPrinciple, radiation conditions and integral formulae for the scatteringof elastic waves”, J. Acoust. Soc. Am. Vol. 59, pp1361-1371 (1976)).Equation (3) can be simplified by identifying proper boundary conditionsfor Σ. In the physical state P, v_(in) and s_(in) are sums of up-goingand down-going waves:v _(in) =v _(in) ^((u)) +v _(in) ^((d))   (4)s _(in) =s _(in) ^((u)) +s _(in) ^((d))   (5)whereas in the hypothetical state Ĥ, {circumflex over (v)}_(in) andŝ_(in) are purely up-going fields:{circumflex over (v)} _(in) ={circumflex over (v)} _(in) ^((u));{circumflex over (v)} _(in) ^((d))=0   (6)ŝ _(in) =ŝ _(in) ^((u)) ;ŝ _(in) ^((d))=0   (7)

These boundary conditions are most conveniently introduced into equation(3) by analysing the problem in the horizontal wavenumber domain, whereup-going and down-going waves and their relation to vertical tractionand particle velocity vectors are explicitly known. Making use ofParceval's identity: $\begin{matrix}{{{\int_{- \infty}^{\infty}\quad{{\mathbb{d}\xi}\quad{f(\xi)}{h(\xi)}}} = {\frac{1}{( {2\pi} )^{2}}{\int_{- \infty}^{\infty}\quad{{\mathbb{d}\kappa}\quad{F( {- \kappa} )}{H(\kappa)}}}}},} & (8)\end{matrix}$where wavenumber vector κ=(k₁, k₂) is conjugate to ξ=(x₁, x₂), equation(3) yields $\begin{matrix}\begin{matrix}{{\overset{\sim}{a}{v_{mn}( {x_{r}❘x_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int{\mathbb{d}{\kappa\lbrack {{{{\hat{S}}_{im}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{V_{in}( {{- \kappa},{z_{r}❘x_{s}}} )}} -} }}}}} \\{ {{{\hat{V}}_{im}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{S_{i\quad n}( {{- \kappa},{z_{r}❘x_{s}}} )}} \rbrack.}\end{matrix} & (9)\end{matrix}$

Going from summation convention to vector notation, equation (9) iswritten $\begin{matrix}\begin{matrix}{{\overset{\sim}{a}{v_{mn}( {x_{r}❘x_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int{\mathbb{d}{\kappa\lbrack {{{{\hat{S}}_{m}^{T}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{V_{n}( {{- \kappa},{z_{r}❘x_{s}}} )}} -} }}}}} \\{ {{{\hat{V}}_{m}^{T}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{S_{n}( {{- \kappa},{z_{r}❘x_{s}}} )}} \rbrack.}\end{matrix} & (10)\end{matrix}$where Ŝ^(T)=(S₁, S₂, S₃) and {circumflex over (V)}^(T)=(V₁, V₂, V₃) arewavenumber domain vertical traction vector and particle velocity vector,respectively, and the superscript ^(T) means transpose.

As shown in the Appendix, since the hypothetical state Ĥ fields Ŝ and{circumflex over (V)} consist of up-going wave modes only, they arerelated asŜ(κ)={circumflex over (L)} _(SV)(κ){circumflex over (V)}(κ),   (11)where {circumflex over (L)}_(SV) is a 3×3 matrix depending on the localmedium parameters ρ, α, β along the receiver spread. The elements of{circumflex over (L)}_(SV) are given in equation (A-19) in the appendix.Inserting equation (11) into equation (10) gives $\begin{matrix}{{{\overset{\sim}{a}{v_{mn}( {x_{r}❘x_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int{{\mathbb{d}\kappa}{{\hat{V}}_{m}^{T}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{R_{n}^{(D)}( {{- \kappa},{z_{r}❘x_{s}}} )}}}}},} & (12)\end{matrix}$whereR ^((D)) ={circumflex over (L)} _(SV) V−S=G ⁻¹ V ^((D))   (13)is interpreted as the overburden response due to scattering in themedium above the receiver level. Furthermore, the vectorV ^((D)) =[V ₁ ^((D)) ,V ₂ ^((D)) ,V ₃ ^((D))]^(T)   (14)contains the elements of the down-going wavemodes on each of theparticle velocity components V₁, V₂, and V₃. Generally, for every shotlocation, V^((D)) is calculated in the slowness or wavenumber domainfrom the particle velocity vector and the vertical traction vectoraccording to equations (A-²⁹), (A-31) and (A-33). The elements of theresulting down-going reflection response are: $\begin{matrix}\begin{matrix}{R_{1}^{(D)} = \frac{2\rho}{p^{2} + {q_{\alpha}q_{\beta}}}} \\{\{ {{\lbrack {q_{\alpha} - {p_{2}^{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}}} \rbrack V_{1}^{(D)}} + {p_{1}p_{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}V_{2}^{(D)}}} \}}\end{matrix} & (15) \\\begin{matrix}{R_{2}^{(D)} = \frac{2\rho}{p^{2} + {q_{\alpha}q_{\beta}}}} \\{\{ {{p_{1}p_{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}V_{1}^{(D)}} + {\lbrack {q_{\alpha} - {p_{1}^{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}}} \rbrack V_{2}^{(D)}}} \}}\end{matrix} & (16) \\{R_{3}^{(D)} = {\frac{2\rho\quad q_{\beta}}{p^{2} + {q_{\alpha}q_{\beta}}}{V_{3}^{(D)}.}}} & (17)\end{matrix}$

Re-introducing the summation convention, we obtain $\begin{matrix}{{\overset{\sim}{a}{v_{mn}( {x_{r}❘x_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int\quad{{\mathbb{d}\kappa}\quad{{\hat{V}}_{im}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{{R_{i\quad n}^{(D)}( {{- \kappa},{z_{r}❘x_{s}}} )}.}}}}} & (18)\end{matrix}$

Making use of Parceval's identity yields in the space domainãv _(mn)(x _(r) | _(s))=∫_(Σ) dS(ξ){circumflex over (v)} _(im)(x|x_(r))r _(in) ^((d))(x|x _(s)),   (19)where r_(in) ^((d)) is the Fourier transform of R_(in) ^((D)).

Equation (19) gives an integral relationship between the field{circumflex over (v)}_(nm) in the hypothetical state Ĥ and the recordedfield v_(nm) in the state P. Note that no information, except location,about the physical source and its radiation characteristics, and noinformation of the properties of the physical overburden above thesurface Σ or the physical subsurface below the receiver level have beenused to derive integral equations (19) for {circumflex over (v)}_(nm).The hypothetical state H fields {tilde over (v)}_(nm) are directlydetermined from reciprocity.

Observe that the inverse of r_(in) ^((d)) in equation (19) can beinterpreted as a multidimensional operator that acts as (i) adeterministic designature operator which removes the physical sourcecharacteristics, and (ii) a deterministic multiple attenuation operatorthat eliminates the overburden response from the physical data.

Eliminating the Incident Wavefield in the Hypothetical State

The desired field {circumflex over (v)}_(nm) in the hypotheticalexperiment can be split into an incident wave field {circumflex over(v)}_(nm) ^((inc)) propagating upwards from the source to the receiver,and the wavefield {circumflex over (v)}_(nm) ^((sc)) scattered upwardsfrom the subsurface: $\begin{matrix}{{\hat{v}}_{n\quad m} = {{\hat{v}}_{n\quad m}^{({inc})} + {{\hat{v}}_{n\quad m}^{({sc})}.}}} & (20)\end{matrix}$

The incident wave field, which propagates in a homogeneous medium andthus is independent of source location, can be shown to be given by:{circumflex over (V)} ^((inc))(z _(r) ⁻ |z _(r))=ãG exp(−iκ·ξ _(r)),  (21)where G is the Green's tensor, equivalent to the tensor derived inequation (A-39) below. Observing that $\begin{matrix}{\begin{matrix}{\begin{matrix}{( {\hat{V}}^{({inc})} )^{T}( {\kappa,{z_{r}^{-}❘x_{r}}} )} \\{R^{(D)}( {{- \kappa},{z_{r}❘x_{s}}} )}\end{matrix} = {\overset{\sim}{a}{G^{T}(\kappa)}{G^{- 1}( {- \kappa} )}{V^{(D)}( {- \kappa} )}{\exp( {{- {\mathbb{i}\kappa}} \cdot \xi_{r}} )}}} \\{= {\overset{\sim}{a}{V^{(D)}( {- \kappa} )}{\exp( {{- {\mathbb{i}\kappa}} \cdot \xi_{r}} )}}}\end{matrix},{and}} & (22) \\{{{v_{mn}^{(d)}( {x_{r}❘x_{s}} )} = {\frac{1}{( {2\pi} )^{2}}{\int{{\mathbb{d}{{\kappa exp}( {{- {\mathbb{i}\kappa}} \cdot \xi_{r}} )}}{V_{mn}^{(D)}( {{- \kappa},{z_{r}❘x_{s}}} )}}}}},} & (23)\end{matrix}$we obtain: $\begin{matrix}{{\overset{\sim}{a}{v_{mn}^{(u)}( {x_{r}❘x_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int{{\mathbb{d}\kappa}\quad{{\hat{V}}_{im}^{({sc})}( {\kappa,{z_{r}^{-}❘x_{r}}} )}{{R_{i\quad n}^{(D)}( {{- \kappa},{z_{r}❘x_{s}}} )}.}}}}} & (24)\end{matrix}$

Making use of Parceval's identity yields in the space domain:ãv _(mn) ^((u))(x _(r) |x _(s))=∫_(Σ) dS(ξ){circumflex over (v)} _(im)^((sc))(x|x _(r))r _(in) ^((D))(x|x _(s)).   (25)

Equation (25) gives the sought-after integral relationship between thescattered field {circumflex over (v)}_(mn) ^((sc)) in the hypotheticalstate Ĥ survey and the up-going and down-going fields v_(mn) ^((u)) andv_(mn) ^((d)) in the state P. The scattered fields {tilde over (v)}_(mn)^((sc)) in the hypothetical state H are determined from the reciprocityrelation.

Equation (25) is a Fredholm integral equation of the first kind for thedesired scattered fields, leading to a system of equations that can besolved for {circumflex over (v)}_(mn) ^((sc)) by keeping the receiverco-ordinate fixed while varying the source co-ordinate.

Wavenumber Domain Solution

Fourier transforming equation (25) over source coordinates ξ, andreceiver coordinates ξ_(r) yields $\begin{matrix}{{\overset{\sim}{a}{V_{mn}^{(U)}( {\kappa_{r},{z_{r}❘\kappa_{s}},z_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int{{\mathbb{d}\kappa}\quad{{\hat{V}}_{im}^{({sc})}( {\kappa,{z_{r}^{-}❘\kappa_{r}},z_{r}} )}{{R_{i\quad n}^{(D)}( {{- \kappa},{z_{r}❘\kappa_{s}},x_{s}} )}.}}}}} & (26)\end{matrix}$

Equation (26) leads to a system of equations that can be solved for{circumflex over (V)}_(mn) ^((sc)) by keeping the wavenumber conjugateto the receiver co-ordinate fixed while varying the wavenumber conjugateto the source co-ordinate. The coupling between positive wavenumbers inthe down-going overburden response field with negative wavenumbers inthe desired field (and vice versa) reflects the autocorrelation processbetween the two fields.

Decomposition into Primary PP, PS, SP, and SS Wave Responses

Following the Betti designature and elastic demultiple step is anelastic wavefield decomposition step that decomposes the multicomponentsource, multicomponent receiver Betti designatured and elasticdemultipled data {circumflex over (v)}^((sc))(x_(r) ^(h)|x_(s) ^(h))given in equation (25) into primary PP, PS, SP, and SS wave responsesthat would be recorded from pressure wave and shear wave sources andreceivers. The receiver and source coordinates for the hypothetical dataare denoted by x_(r) ^(h) and x_(s) ^(h), respectively, where x_(r) ^(h)is a point on the surface Σ and x_(s) ^(h)=x_(r) is infinitesimallybelow. Let matrix Û contain up-going potentials at the receiver fromsources of force in direction i, $\begin{matrix}{\hat{U} = {\begin{pmatrix}{\hat{U}}_{P1} & {\hat{U}}_{P2} & {\hat{U}}_{P3} \\{\hat{U}}_{SV1} & {\hat{U}}_{SV2} & {\hat{U}}_{SV3} \\{\hat{U}}_{SH1} & {\hat{U}}_{SH2} & {\hat{U}}_{SH3}\end{pmatrix}.}} & (27)\end{matrix}$

Let U be the matrix containing up-going P, SV, and SH wave potentials atthe receiver excited from P, SV, and SH wave sources, $\begin{matrix}{U = {\begin{pmatrix}U_{PP} & U_{PSV} & U_{PSH} \\U_{SVP} & U_{SVSV} & U_{SVSH} \\U_{SHP} & U_{SHSV} & U_{SHSH}\end{pmatrix}.}} & (28)\end{matrix}$

This elastic source-receiver wavefield decomposition may now be dividedinto two computational operations. First, the multicomponent data{circumflex over (v)}^((sc)) are Fourier transformed over horizontalreceiver co-ordinates ξ_(r) ^(h). In the wavenumber domain κ_(r),receiver side wavefield decomposition into upgoing pressure and shearwavefields is achieved by left-multiplying {circumflex over (V)}^((sc))by a receiver decomposition matrix R,Û(κ_(r) , z _(r) ^(h) |x _(s) ^(h))=R(κ_(r)){circumflex over (V)}^((sc))(κ_(r) , z _(r) ^(h) |x _(s) ^(h)),   (29)where $\begin{matrix}{{R = {\frac{1}{p^{2} + {q_{\alpha}q_{\beta}}}\begin{pmatrix}p_{1} & p_{2} & {- q_{\beta}} \\{- \frac{p_{1}q_{\alpha}}{p}} & {- \frac{p_{2}q_{\alpha}}{p}} & {- p} \\{{- \frac{p_{2}}{p}}( {p^{2} + {q_{\alpha}q_{\beta}}} )} & {\frac{p_{1}}{p}( {p^{2} + {q_{\alpha}q_{\beta}}} )} & 0\end{pmatrix}}},} & (30)\end{matrix}$with p=κ_(r)/ω. An inverse Fourier transform gives upgoing pressure andshear wavefields measured at the receiver Û(x_(r) ^(h)|x_(s) ^(h)).Second, the multicomponent data Û(x_(r) ^(h)|x_(s) ^(h)) are Fouriertransformed over horizontal source co-ordinates ξ_(s) ^(h). In thewavenumber domain κ_(s), source side wavefield decomposition intodown-going pressure and shear source wavefields is achieved byright-multiplying the data by a source decomposition operator S,U(x _(r) ^(h) |κ _(s) ,z _(s) ^(h))=Û(x _(r) ^(h) |κ _(s) ,z _(s)^(h))S(κ_(s)),   (31)where $\begin{matrix}{{S = {\frac{{{- i}\quad\omega^{- 1}\rho}\quad}{p^{2} + {q_{\alpha}q_{\beta}}}\begin{pmatrix}p_{1} & {- \frac{p_{1}q_{\alpha}}{p}} & {{- \frac{p_{2}q_{\alpha}}{p}}{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}} \\p_{2} & {- \frac{p_{2}q_{\alpha}}{p}} & {\frac{p_{1}q_{\alpha}}{p}{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}} \\{- q_{\beta}} & {- p} & 0\end{pmatrix}}},} & (32)\end{matrix}$with p=κ_(s)/ω. An inverse Fourier transform gives data U equivalent todata from a hypothetical survey with overburden absent, with singlecomponent pressure and shear wave sources, and single component pressureand shear wave receivers.

In (κ_(r), κ_(s)) space, the source-receiver decomposition isU(κ_(r) ,z _(r) ^(h) |κ _(s) ,z _(s) ^(h))=R(κ_(r)){circumflex over (V)}^((sc))(κ_(r) ,z _(r) ^(h) |κ _(s) ,z _(s) ^(h))S(κ_(s)),   (33)

Observe that the decomposition operators are simply matrices which aremultiplied by the data. The required parameters for the decompositionare the density, and P-wave and S-wave velocities at the receiver depthlevel.

Betti Deconvolution and Source-Receiver Wavefield Decomposition:Horizontally Layered Medium (“1.5D Medium”)

In a horizontally layered medium, the response is laterallyshift-invariant with respect to horizontal source location ξ_(s), andthe expressions therefore simplify significantly. For this reason, ahorizontally layered medium ls sometimes referred to as a “1.5D medium”.Assuming a horizontally layered medium often yields satisfactory resultsfor real data and realistic Earth structures.

Consider equation (24), where we may set ξ_(s)=0 and{circumflex over (V)} _(im) ^((sc))(κ,z _(r) ⁻ |x _(r))={circumflex over(V)} _(im) ^((sc))(κ,z _(r) ⁻|ξ=0,z _(r))exp(−iκ·ξ _(r)).   (34)

By Fourier transforming equation (24) with respect to ξ_(r) andinterchanging integrals, we find $\begin{matrix}{{\overset{\sim}{a}{V_{mn}^{U}( {\kappa_{r},{{z_{r}❘\xi_{s}} = 0},z_{s}} )}} = {\frac{1}{( {2\pi} )^{2}}{\int_{- \infty}^{\infty}\quad{{\mathbb{d}\kappa}\quad{{\hat{V}}_{im}^{({sc})}( {\kappa,{{z_{r}❘\xi} = 0},z_{r}} )} \times {R_{in}^{(D)}( {{- \kappa},{{z_{r}^{-}❘\xi_{s}} = 0},z_{s}} )}{\int_{- \infty}^{\infty}\quad{{\mathbb{d}\xi_{r}}{{\exp( {{- i}\quad{\xi_{r} \cdot ( {\kappa + \kappa_{r}} )}} )}.}}}}}}} & (35)\end{matrix}$

The last integral is recognised as the Dirac delta function. Performingthe integration over wavenumbers, using the Dirac delta functionproperty, and renaming κ_(r) by κ, we obtainãV _(mn) ^(U)(κ,z _(r)|ξ_(s)=0,z _(s))={circumflex over (V)} _(im)^((sc))(−κ,z _(r) ⁻|ξ=0,z _(r))R _(in) ^((D))(κ,z _(r) ⁻|ξ_(s)=0,z _(s))  (36)

Writing this result in terms of matrices, yields $\begin{matrix}{\begin{pmatrix}{\hat{V}}_{11}^{({sc})} & {\hat{V}}_{21}^{({sc})} & {\hat{V}}_{31}^{({sc})} \\{\hat{V}}_{12}^{({sc})} & {\hat{V}}_{22}^{({sc})} & {\hat{V}}_{32}^{({sc})} \\{\hat{V}}_{13}^{({sc})} & {\hat{V}}_{23}^{({sc})} & {\hat{V}}_{33}^{({sc})}\end{pmatrix}_{({{- \kappa},{z_{r}^{-}❘z_{r}}})} = {{\overset{\sim}{a}\begin{pmatrix}V_{11}^{(U)} & V_{12}^{(U)} & V_{13}^{(U)} \\V_{21}^{(U)} & V_{22}^{(U)} & V_{23}^{(U)} \\V_{31}^{(U)} & V_{32}^{(U)} & V_{33}^{(U)}\end{pmatrix}}_{({\kappa,{z_{r}❘z_{r}}})}{\begin{pmatrix}R_{11}^{(D)} & R_{12}^{(D)} & R_{13}^{(D)} \\R_{21}^{(D)} & R_{22}^{(D)} & R_{23}^{(D)} \\R_{31}^{(D)} & R_{32}^{(D)} & R_{33}^{(D)}\end{pmatrix}^{- 1}}_{({\kappa,{z_{r}❘z_{r}}})}}} & (37)\end{matrix}$

Using the field properties $\begin{matrix}{{{\hat{V}}_{ij}^{({sc})}( {- \kappa} )} = \{ {\begin{matrix}{{\hat{V}}_{ij}^{({sc})}(\kappa)} & {i = j} \\{- {{\hat{V}}_{ij}^{({sc})}(\kappa)}} & {i \neq j}\end{matrix}.} } & (38)\end{matrix}$it is evident that the scattered part of the desired field is obtainedby generalised spectral deconvolution: $\begin{matrix}{\begin{pmatrix}{\hat{V}}_{11}^{({sc})} & {- {\hat{V}}_{21}^{({sc})}} & {- {\hat{V}}_{31}^{({sc})}} \\{- {\hat{V}}_{12}^{({sc})}} & {\hat{V}}_{22}^{({sc})} & {- {\hat{V}}_{32}^{({sc})}} \\{- {\hat{V}}_{13}^{({sc})}} & {- {\hat{V}}_{23}^{({sc})}} & {\hat{V}}_{33}^{({sc})}\end{pmatrix}_{({\kappa,{z_{r}^{-}❘z_{r}}})} = {{\overset{\sim}{a}\begin{pmatrix}V_{11}^{(U)} & V_{12}^{(U)} & V_{13}^{(U)} \\V_{21}^{(U)} & V_{22}^{(U)} & V_{23}^{(U)} \\V_{31}^{(U)} & V_{32}^{(U)} & V_{33}^{(U)}\end{pmatrix}}_{({\kappa,{z_{r}❘z_{r}}})}{\begin{pmatrix}R_{11}^{(D)} & R_{12}^{(D)} & R_{13}^{(D)} \\R_{21}^{(D)} & R_{22}^{(D)} & R_{23}^{(D)} \\R_{31}^{(D)} & R_{32}^{(D)} & R_{33}^{(D)}\end{pmatrix}^{- 1}}_{({\kappa,{z_{r}❘z_{r}}})}}} & (39)\end{matrix}$

It follows from equation (39) that the components of the desired fieldare obtained by deterministic spectral deconvolution between the fielditself and the reflection response of the overburden containing thedown-going part of the particle velocity vector.

The next step is to decompose the Betti deconvolved data into primaryPP, PS, SP, and SS wave responses. This is carried out as described inthe section “Decomposition into primary PP, PS, SP, and SS waveresponses”.

Observe that the 1.5D Betti deconvolution scheme followed bysource-receiver decomposition may be implemented as τ−p orfrequency-wavenumber domain algorithms. In the τ−p domain, a jointdesignature, multiple attenuation and source-receiver decompositionprocess is performed for each p-trace. In the frequency-wavenumberdomain, the process is performed for each combination of frequency andwavenumber.

Applications for the Present Invention

As noted above, the present invention can be applied to a number ofdifferent seismic survey configurations. Three survey scenarios ofparticular interest will be discussed briefly.

1. Land Surface Seismic Acquisition

The invention may be applied to a 3C×3C (or 9C) land seismic surveywhere three orthogonal source motions (two horizontal and one vertical)are generated separately, and where the seismic wavefield is recorded bythree orthogonal geophones, measuring the two horizontal and thevertical components of the particle velocity vector.

For measurements located at the free surface, the up/down separationsimplifies considerably: all tractions S₁, S₂, and S₃ vanish so that thetraction terms vanish in equations (A-28) to (A-33) below.

2. Seabed Seismic Acquisition

The invention may be applied to seismic data acquired in a 3C×4C (or12C) seabed seismic experiment where three orthogonal source motions(two horizontal and one vertical) are generated separately, either onthe ocean bottom or in the water column, and where the seismic wavefieldis recorded at a receiver location by three orthogonal geophonesdeployed just below the sea floor and one hydrophone deployed just abovethe sea floor.

For measurements located at the fluid/solid interface of the seabed, theup/down separation simplifies considerably: tractions S₁ and S₂ vanishso that the corresponding traction terms vanish in equations (A-28) to(A-33) below. The hydrophone measurement on the seafloor is the negativeof S₃.

It is worthwhile to note that sources that are equivalent to sources offorce can be generated in the water by using, for instance, aconventional airgun source. For instance, the response corresponding toa vertical point-force can be obtained by acquiring two records with thesame source at the same location but at slightly different depths. Afinite-difference approximation then allows us to compute the responsedue to a pressure gradient source in the vertical direction. Theequation of motion (Newton's second law) shows that this is equivalentto a vertical point force. The same is true for the two horizontaldirections, as well.

3. Borehole Seismic Acquisition

The 3C×6C (or 18C) borehole seismic experiment where the threeorthogonal source motions are generated separately, either on land, onthe sea floor, or in the water column, for which the three components ofparticle velocity and the three components of vertical traction areknown.

In this surveying arrangement, there are no boundary conditions thatsimplify the up/down separation step. However, by adding hydrophones inthe receiver package, additional constraints can be obtained (ahydrophone measures divergence of particle motion).

Wavenumber Domain Basic Relationships Between the Particle Velocity andVertical Traction Vectors and Upgoing and Downgoing Wave Vectors

This section considers a horizontally layered elastic earth. In asource-free region elastic wave propagation is described by the equationof motion and Hooke's law (the elastic constitutive relation). These canbe written as a system of first-order ordinary differential equations ofthe form∂₃B=−iωAB   (A-1)where the field vector B is definedB=(V ^(T) ,S ^(T))^(T)   (A-2)with particle velocity vector V^(T)=[V₁, V₂, V₃] and vertical tractionvector S^(T)=[S₁, S₂, S₃]. The elastic system matrix A has the form$\begin{matrix}{A = \begin{pmatrix}0 & 0 & p_{1} & \frac{1}{\mu} & 0 & 0 \\0 & 0 & p_{2} & 0 & \frac{1}{\mu} & 0 \\{\frac{\lambda}{\lambda + {2\mu}}p_{1}} & {\frac{\lambda}{\lambda + {2\mu}}p_{2}} & 0 & 0 & 0 & \frac{1}{\lambda + {2\mu}} \\{\rho - {\theta\quad p_{1}^{2}} - {\mu\quad p_{\sigma}p_{\sigma}}} & {{- \theta}\quad p_{1}p_{2}} & 0 & 0 & 0 & {\frac{\lambda}{\lambda + {2\mu}}p_{1}} \\{{- \theta}\quad p_{1}p_{2}} & {\rho - {\theta\quad p_{2}^{2}} - {\mu\quad p_{\sigma}p_{\sigma}}} & 0 & 0 & 0 & {\frac{\lambda}{\lambda + {2\mu}}p_{2}} \\0 & 0 & \rho & p_{1} & p_{2} & 0\end{pmatrix}} & ( {A\text{-}3} )\end{matrix}$where θ=μ(3λ+2 μ)/(μ+2 μ) and p²=p₁ ²+p₂ ².

For notational convenience, the explicit dependence of differentquantities on frequency, wavenumber, depth, etc., is omitted. Forinstance, the particle velocity vector recorded at depth x₃,v(ξ,x₃,ω,x_(s)), due to a point source at location x_(s) is in thewavenumber domain denoted V or V(x₃) with the understandingV=V(x₃)==V(κ,x₃,ω,x_(s)). When required we will show dependency onhorizontal slowness vector p=κ/ω=(p₁,p₂)

Up-Going and Down-Going Waves

For the docomposition of the elastic field into up- and down-going wavesin a layered earth it is necessary to find the eigenvalues andeigenvectors of the system matrix A for given wavenumbers andfrequencies. The field vector B can be decomposed into up-(U) anddown-going (D) wavesW=[U ^(T) ,D ^(T)]^(T)   (A-4)where U^(T)=[U_(P), U_(SV), U_(SH)] and D^(T)=[D_(P), D_(SV), D_(SH)],by the linear transformationB=LW   (A-5)where L is the local eigenvector matrix of A (i.e., each column of L isan eigenvector). Equation (A-5) describes composition of the wavefield Bfrom its upgoing and downgoing constituents.

Given the inverse eigenvector matrix L⁻¹, the up- and downgoing wavescan be computed by evaluatingW=L⁻¹B.   (A-6)

Equation (A-6) describes decomposition of the wavefield B into upgoingand downgoing P- and S-waves. After some straightforward, but tediouscalculations, the composition matrix is obtained $\begin{matrix}{{L = {{L(p)} = \begin{pmatrix}{L_{VU}(p)} & {- {L_{VU}( {- p} )}} \\{L_{SU}(p)} & {L_{SU}( {- p} )}\end{pmatrix}}},} & ( {A\text{-}7} )\end{matrix}$and the decomposition matrix $\begin{matrix}{{L^{- 1} = {{L^{- 1}(p)} = \begin{pmatrix}{L_{SU}^{T}(p)} & {L_{VU}^{T}(p)} \\{- {L_{SU}^{T}( {- p} )}} & {L_{VU}^{T}( {- p} )}\end{pmatrix}}},{where}} & ( {A\text{-}8} ) \\{{L_{VU} = {\frac{1}{\sqrt{2}}\begin{pmatrix}{{- p_{1}}\frac{1}{\sqrt{\rho\quad q_{\alpha}}}} & {\frac{p_{1}}{p}\sqrt{\frac{q_{\beta}}{\rho}}} & {\frac{p_{2}}{p}\frac{1}{\sqrt{\mu\quad q_{\beta}}}} \\{{- p_{2}}\frac{1}{\sqrt{\rho\quad q_{\alpha}}}} & {\frac{p_{2}}{p}\sqrt{\frac{q_{\beta}}{\rho}}} & {{- \frac{p_{1}}{p}}\frac{1}{\sqrt{\mu\quad q_{\beta}}}} \\\sqrt{\frac{q_{\alpha}}{\rho}} & {p\frac{1}{\sqrt{\rho\quad q_{\beta}}}} & 0\end{pmatrix}}},} & ( {A\text{-}9} ) \\{L_{SU} = {\frac{1}{\sqrt{2}}{\begin{pmatrix}{{- 2}\mu\quad p_{1}\sqrt{\frac{q_{\alpha}}{\rho}}} & {\frac{p_{1}}{p}( {\rho - {2\mu\quad p^{2}}} )\frac{1}{\sqrt{\rho\quad q_{\beta}}}} & {\frac{p_{2}}{p}\sqrt{\mu\quad q_{\beta}}} \\{{- 2}\mu\quad p_{2}\sqrt{\frac{q_{\alpha}}{\rho}}} & {\frac{p_{2}}{p}( {\rho - {2\mu\quad p^{2}}} )\frac{1}{\sqrt{\rho\quad q_{\beta}}}} & {{- \frac{p_{1}}{p}}\sqrt{\mu\quad q_{\beta}}} \\{( {\rho - {2\mu\quad p^{2}}} )\frac{1}{\sqrt{\rho\quad q_{\alpha}}}} & {2\mu\quad p\sqrt{\frac{q_{\beta}}{\rho}}} & 0\end{pmatrix}.}}} & ( {A\text{-}10} )\end{matrix}$

In a source-free homogeneous solid upgoing and downgoing waves satisfythe differential equations:∂₃ U _(P) =−iωq _(a) U _(P),   (A-11)∂₃ U _(SV) =−iωq _(β) U _(SV),   (A-12)∂₃ U _(SH) =−iωq _(β) U _(SH),   (A-13)∂₃D_(P)=iωq_(a)D_(P),   (A-14)∂₃D_(SV)=iωq_(β)D_(SV),   (A-15)∂₃D_(SH)=iωq_(β)D_(SH).   (A-16)Traction-Particle Velocity Vector Relationship When Down-Going WavesVanishUsing the relation (A-5) withD=0one obtains a simple relationship between the vertical traction vectorand the particle velocity vector,S={circumflex over (L)}_(SV)V,   (A-17)where{circumflex over (L)}_(SV)={circumflex over (L)}_(SU){circumflex over(L)}_(VU) ⁻¹.   (A-18)$\begin{matrix}{{\hat{L}}_{SV} = {\frac{\rho}{p^{2} + {q_{\alpha}q_{\beta}}}\begin{pmatrix}{q_{\alpha} - {p_{2}^{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}}} & {p_{1}p_{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}} & {p_{1}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack} \\{p_{1}p_{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}} & {q_{\alpha} - {p_{1}^{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}}} & {p_{2}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack} \\{- {p_{1}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}} & {- {p_{2}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}} & q_{\beta}\end{pmatrix}}} & ( {A\text{-}19} )\end{matrix}$Upgoing and Downgoing Waves Scaled as Particle Velocity and VerticalTraction

Upgoing and downgoing P, SV, and SH waves are not uniquely defined. Wemay scale the components so that they have dimensions of particlevelocity or traction in such a way that their sum gives a component ofparticle velocity or a component of vertical traction. Let V_(i) ^((U)^(P) ⁾ denote the upgoing P-waves on V_(i), V_(i) ^((U) ^(SV) ⁾ denotethe upgoing SV-waves on V_(i), etc. Then,V _(i) =V _(i) ^((U) ^(P) ⁾ +V _(i) ^((U) ^(SV) ⁾ +V _(i) ^((U) ^(SH) ⁾+V _(i) ^((D) ^(P) ⁾ +V _(i) ^((D) ^(SV) ⁾ +V _(i) ^((D) ^(SH) ⁾with similar equation for vertical traction. Furthermore, it is possibleto sum the upgoing components into a total upgoing component, andlikewise for the downgoing components. V_(i) ^((U)) is defined as thesum of upgoing waves on V_(i),V _(i) ^((U)) =V _(i) ^((U) ^(P) ⁾ +V _(i) ^((U) ^(SV) ⁾ +V _(i) ^((U)^(SH) ⁾,and V_(i) ^((D)) as the sum of downgoing waves on V_(i),V _(i) ^((D)) =V _(i) ^((D) ^(P) ⁾ +V _(i) ^((D) ^(SV) ⁾ +V _(i) ^((D)^(SH) ⁾.

The total upgoing waves on S_(i) areS _(i) ^((U)) =S _(i) ^((U) ^(P) ⁾ +S _(i) ^((U) ^(SV) ⁾ +S _(i) ^((U)^(SH) ⁾,and the total downgoing waves on S_(i),S _(i) ^((D)) =S _(i) ^((D) ^(P) ⁾ +S _(i) ^((D) ^(SV) ⁾ +S _(i) ^((D)^(SH) ⁾.

Further, defining the vectorsV ^((U)) =[V ₁ ^((U)) ,V ₂ ^((U)) ,V ₃ ^((U))]^(T)   (A-20)V ^((D)) =[V ₁ ^((D)) ,V ₂ ^((D)) ,V ₃ ^((D))]^(T)   (A-21)S ^((U)) =[S ₁ ^((U)) ,S ₂ ^((U)) ,S ₃ ^((U))]^(T)   (A-22)S ^((D)) =[S ₁ ^((D)) ,S ₂ ^((D)) ,S ₃ ^((D))]^(T)   (A-23)it follows that these vectors are related to the originally definedvectors of upgoing and downgoing waves, U and D, respectively, asV ^((U))(p)=L _(VU)(p)U(p)   (A-24)V ^((D))(p)=−L _(VU)(−p)D(p)   (A-25)S ^((U))(p)=L _(SU)(p)U(p)   (A-26)S ^((D))(p)=L _(SU)(−p)D(p)   (A-27)V^((D)) and V^((U)) enter the integral relationship between the physicalexperiment and the hypothetical experiment as known field vectors.

Explicitly, the total upgoing and downgoing waves on V₁ are given fromV₁, V₃, and S₁ as $\begin{matrix}{{V_{1}^{(U)} = {{\frac{1}{2}V_{1}} - {{\frac{p_{1}}{2q_{\alpha}}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}V_{3}} + {\frac{1}{2\rho\quad q_{\alpha}}( {p^{2} + {q_{\alpha}q_{\beta}}} )S_{1}}}},} & ( {A\text{-}28} ) \\{V_{1}^{(D)} = {{\frac{1}{2}V_{1}} + {{\frac{p_{1}}{2q_{\alpha}}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}V_{3}} - {\frac{1}{2\rho\quad q_{\alpha}}( {p^{2} + {q_{\alpha}q_{\beta}}} ){S_{1}.}}}} & ( {A\text{-}29} )\end{matrix}$

Note that total up-going and down-going waves on V₁ do not depend on V₂,S₂, and S₃. Further, total up-going and down-going waves on V₂ are givenfrom V₂, V₃, and S₁ as $\begin{matrix}{{V_{2}^{(U)} = {{\frac{1}{2}V_{2}} - {{\frac{p_{2}}{2q_{\alpha}}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}V_{3}} + {\frac{1}{2\rho\quad q_{\alpha}}( {p^{2} + {q_{\alpha}q_{\beta}}} )S_{1}}}},} & ( {A\text{-}30} ) \\{{V_{2}^{(D)} = {{\frac{1}{2}V_{2}} + {{\frac{p_{2}}{2q_{\alpha}}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}V_{3}} - {\frac{1}{2\rho\quad q_{\alpha}}( {p^{2} + {q_{\alpha}q_{\beta}}} )S_{1}}}},} & ( {A\text{-}31} )\end{matrix}$

and are thus independent of V₁, S₂, and S₃. Total up-going anddown-going waves on V₃ are given from V₁, V₂, V₃, and S₃ as$\begin{matrix}{{V_{3}^{(U)} = {{\frac{1}{2}V_{3}} + {\frac{1}{2\rho\quad q_{\beta}}( {p^{2} + {q_{\alpha}q_{\beta}}} )S_{3}} + {{\frac{1}{2q_{\beta}}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}( {{p_{1}V_{1}} + {p_{2}V_{2}}} )}}},} & ( {A\text{-}32} ) \\{V_{3}^{(D)} = {{\frac{1}{2}V_{3}} - {\frac{1}{2\rho\quad q_{\beta}}( {p^{2} + {q_{\alpha}q_{\beta}}} )S_{3}} - {{\frac{1}{2q_{\beta}}\lbrack {1 - {2{\beta^{2}( {p^{2} + {q_{\alpha}q_{\beta}}} )}}} \rbrack}{( {{p_{1}V_{1}} + {p_{2}V_{2}}} ).}}}} & ( {A\text{-}33} )\end{matrix}$Derivation of Matrix Relations

Equation (A-18) and symmetry relations of (A-19) yields $\begin{matrix}{\begin{matrix}{{\lbrack {{{\hat{L}}_{SV}V} - S} \rbrack(p)} = {{{L_{VU}^{- T}( {- p} )}{L_{SU}^{T}( {- p} )}{V(p)}} - {S(p)}}} \\{= {{L_{VU}^{- T}( {- p} )}\lbrack {{{L_{SU}^{T}( {- p} )}{V(p)}} - {{L_{VU}^{T}( {- p} )}{S(p)}}} \rbrack}}\end{matrix}.} & ( {A\text{-}34} )\end{matrix}$

The last term in the square brackets is identified by use of equations(A-6) and (A-8) as the downgoing wave vector calculated from V and Saccording toD(p)=L _(SU) ^(T)(−p)V(p)−L _(VU) ^(T)(−p)S(p)   (A-35)

In the previous subsection, it was shown that D was related to thedown-going particle ave vector V(D) by equations (A-28) and (A-29).Inverting this equation givesD(p)=−L _(VU) ⁻¹(−p)V ^((D))(p),   (A-36)which is substituted into equation (A-34) to give $\begin{matrix}{\begin{matrix}{{\lfloor {{{\hat{L}}_{SV}V} - S} \rfloor(p)} = {{L_{VU}^{- T}( {- p} )}{L_{VU}^{- 1}( {- p} )}{V^{(D)}(p)}}} \\{{= {\lbrack {{L_{VU}( {- p} )}{L_{VU}^{T}( {- p} )}} \rbrack^{1}{V^{(D)}(p)}}},} \\{= {{G^{- 1}(p)}{V^{(D)}(p)}}}\end{matrix}{where}} & ( {A\text{-}37} ) \\{G^{- 1} = {\frac{1}{2}{( {L_{VU}L_{VU}^{T}} )^{- 1}.}}} & ( {A\text{-}38} )\end{matrix}$

Explicitly, it is found that $\begin{matrix}{G^{- 1} = {\frac{2\rho}{p^{2} + {q_{\alpha}q_{\beta}}}{\begin{pmatrix}{q_{\alpha} - {p_{2}^{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}}} & {p_{1}p_{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}} & 0 \\{p_{1}p_{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}} & {q_{\alpha} - {p_{1}^{2}{\beta^{2}( {q_{\alpha} - q_{\beta}} )}}} & 0 \\0 & 0 & q_{\beta}\end{pmatrix}.}}} & ( {A\text{-}39} )\end{matrix}$

In the theory given above it has been assumed that the sources in theidealised survey are point sources of force. However, the invention isnot limited to the use of sources in the idealised survey which arepoint sources of force, although the detailed theory given aboverequires slight modification to accommodate sources in the idealisedsurvey that are not point sources of force.

In the detailed description given above the designature/demultipleoperator is derived using Betti's theorem. The invention is not howeverlimited to deriving the designature/demultiple operator using Betti'stheorem. The designature/demultiple operator of the invention may bederived from the elastodynamic wave equation or a representation of theelastodynamic wave equation. Betti's theorem is an integralrepresentation of the elastodynamic wave equation, and examples of otherrepresentations of the elastodynamic wave equation include the elasticKirchhoff integral or the elastodynamic representation theorem, and theinverse scattering series method.

FIG. 12 is a schematic block diagram of an apparatus 15 according to thepresent invention. The apparatus is able to carry out a method accordingto the present invention.

The apparatus 15 comprises a programmable data processor 16 with aprogram memory 17, for instance in the form of a read only memory ROM,storing a program for controlling the data processor 17 to perform amethod of the invention. The system further comprises non-volatileread/write memory 18 for storing, for example, any data which must beretained in the absence of power supply. A “working” or “scratchpad”memory for the data processor is provided by a random access memory(RAM) 19. An input device 20 is provided, for instance for receivinguser commands and data. An output device 21 is provided, for instancefor displaying information relating to the progress and result of themethod. The output device may be, for example, a printer, a visualdisplay unit or an output memory.

Seismic data for processing according to a method of the invention maybe supplied via the input device 20 or may optionally be provided by amachine-readable store 22.

The program for operating the system and for performing the methoddescribed hereinbefore is stored in the program memory 17, which may beembodied as a semi-conductor memory, for instance of the well-known ROMtype. However, the program may be stored in any other suitable storagemedium, such as magnetic data carrier 17 a (such as a “floppy disc”) orCD-ROM 17 b.

1. A method of processing multi-component seismic data, the data havingbeen acquired by emitting multi-component seismic energy at a sourcelocation; and acquiring seismic data at a multi-component seismicreceiver located at a greater depth than the source location, the methodcomprising the steps of: decomposing the acquired seismic data intoup-going constituents and down-going constituents; and calculating ade-signature and de-multiple operator from the down-going constituentsof the acquired seismic data and from properties of the mediumsurrounding the receiver.
 2. A method as claimed in claim 1 andcomprising the further step of processing the acquired seismic datausing the de-signature and de-multiple operator thereby to attenuate orremove seismic events arising from multiple reflections.
 3. A method asclaimed in claim 2 wherein the step of processing the acquired seismicdata comprises processing the up-going constituents of the acquiredseismic data using the de-signature and de-multiple operator.
 4. Amethod as claimed in claim 2 wherein the step of processing the acquiredseismic data further comprises selecting a desired seismic signature forthe source.
 5. A method as claimed in claim 2 and comprising the furtherstep of decomposing the processed seismic data into P-wave and/or S-wavedata.
 6. A method as claimed in claim 1 and comprising the step ofdecomposing the acquired seismic data into P-wave and/or S-wave data,the step of decomposing the acquired seismic data into P-wave and/orS-wave data being performed before the step of decomposing the seismicdata into up-going constituents and down-going constituents.
 7. A methodas claimed in claim 5 wherein the step of decomposing the seismic datainto P-wave and/or S-wave data comprises decomposing the data at thereceiver side.
 8. A method as claimed in claim 5 wherein the step ofdecomposing the seismic data into P-wave and/or S-wave data comprisesdecomposing the data at the source side.
 9. A method of processingmulti-component seismic data, the data having been acquired by emittingmulti-component seismic energy at a source location; and acquiringseismic data at a multi-component seismic receiver located at a greaterdepth than the source location, the method comprising the steps of:decomposing the acquired seismic data into up-going constituents anddown-going constituents; and calculating a de-signature operator fromthe initial down-going constituents of the acquired seismic data andfrom properties of the medium surrounding the receiver.
 10. A method ofprocessing multi-component seismic data, the data having been acquiredby emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location, the method comprising thesteps of: decomposing the acquired seismic data into PP-constituents,PS-constituent, SP-constituents and SS-constituents; and calculating ade-signature and de-multiple operator from the constituents of theacquired seismic data and from properties of the medium surrounding thereceiver.
 11. A method as claimed in claim 10 and comprising the furtherstep of processing the acquired seismic data using the de-signature andde-multiple operator thereby to attenuate or remove seismic eventsarising from multiple reflections.
 12. A method of calculating ade-signature and de-multiple operator for multi-component seismic data,the data having been acquired by emitting multi-component seismic energyat a source location; and acquiring seismic data at a multi-componentseismic receiver located at a greater depth than the source location,the method comprising the steps of: decomposing the acquired seismicdata into up-going constituents and down-going constituents; andcalculating a de-signature and de-multiple operator from the down-goingconstituents of the acquired seismic data and from properties of themedium surrounding the receiver.
 13. A method of calculating ade-signature operator for multi-component seismic data, the data havingbeen acquired by emitting multi-component seismic energy at a sourcelocation; and acquiring seismic data at a multi-component seismicreceiver located at a greater depth than the source location, the methodcomprising the steps of: decomposing the acquired seismic data intoup-going constituents and down-going constituents; and calculating ade-signature operator from the initial down-going constituents of theacquired seismic data and from properties of the medium surrounding thereceiver.
 14. A method of calculating a de-signature and de-multipleoperator for multi-component seismic data, the data having been acquiredby emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location, the method comprising thesteps of: decomposing the acquired seismic data into PP-constituents,PS-constituent, SP-constituents and SS-constituents; and calculating ade-signature and de-multiple operator from the constituents of theacquired seismic data and from properties of the medium surrounding thereceiver.
 15. A method as claimed in claim 1 and further comprising thesteps of: emitting multi-component seismic energy at a source location;and acquiring seismic data at a multi-component seismic receiver locatedat a greater depth than the source location.
 16. An apparatus forprocessing multi-component seismic data, the date having been acquiredby emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location, the apparatus comprising:means for decomposing the acquired seismic data into up-goingconstituents and down-going constituents; and means for calculating ade-signature and de-multiple operator from the down-going constituentsof the acquired seismic data and from properties of the mediumsurrounding the receiver.
 17. An apparatus as claimed in claim 16 andfurther comprising means for processing the acquired seismic data usingthe de-signature and de-multiple operator thereby to attenuate or removeseismic events arising from multiple reflections.
 18. An apparatus asclaimed in claim 17 and further comprising means for processing theup-going constituents of the acquired seismic data using thede-signature and de-multiple operator.
 19. An apparatus as claimed inclaim 16 and further comprising means for selecting a desired seismicsignature for the source.
 20. An apparatus as claimed in claim 17 andfurther comprising means for decomposing the processed seismic data intoP-wave and/or S-wave data.
 21. An apparatus as claimed in claim 20 andfurther comprising means for decomposing the seismic data into P-waveand/or S-wave data at the receiver side.
 22. An apparatus as claimed inclaim 20 and further comprising means for decomposing the seismic datainto P-wave and/or S-wave data at the source side.
 23. An apparatus forprocessing multi-component seismic data, the data having been acquiredby emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location, the apparatus comprising:means for decomposing the acquired seismic data into up-goingconstituents and down-going constituents; and means for calculating ade-signature operator from the initial down-going constituents of theacquired seismic data and from properties of the medium surrounding thereceiver.
 24. An apparatus for processing multi-component seismic data,the data having been acquired by emitting multi-component seismic energyat a source location; and acquiring seismic data at a multi-componentseismic receiver located at a greater depth than the source location,the apparatus comprising: means for decomposing the acquired seismicdata into PP-constituents, PS-constituent, SP-constituents andSS-constituents; and means for calculating a de-signature andde-multiple operator from the constituents of the acquired seismic dataand from properties of the medium surrounding the receiver.
 25. A methodas claimed in claim 9 and comprising the further step of processing theacquired seismic data using the de-signature and de-multiple operatorthereby to attenuate or remove seismic events arising from multiplereflections.
 26. An apparatus for calculating a de-signature andde-multiple operator for multi-component seismic data, the data havingbeen acquired by emitting multi-component seismic energy at a sourcelocation; and acquiring seismic data at a multi-component seismicreceiver located at a greater depth than the source location, theapparatus comprising: means for decomposing the acquired seismic datainto up-going constituents and down-going constituents; and means forcalculating a de-signature and de-multiple operator from the down-goingconstituents of the acquired seismic data and from properties of themedium surrounding the receiver.
 27. An apparatus for calculating ade-signature operator for multi-component seismic data, the data havingbeen acquired by emitting multi-component seismic energy at a sourcelocation; and acquiring seismic data at a multi-component seismicreceiver located at a greater depth than the source location, theapparatus comprising: means for decomposing the acquired seismic datainto up-going constituents and down-going constituents; and means forcalculating a de-signature operator from the initial down-goingconstituents of the acquired seismic data and from properties of themedium surrounding the receiver.
 28. An apparatus for calculating ade-signature and de-multiple operator for multi-component seismic data,the data having been acquired by emitting multi-component seismic energyat a source location; and acquiring seismic data at a multi-componentseismic receiver located at a greater depth than the source location,the apparatus comprising: means for decomposing the acquired seismicdata into PP-constituents, PS-constituent, SF-constituents andSS-constituents; and means for calculating a de-signature andde-multiple operator from the constituents of the acquired seismic dataand from properties of the medium surrounding the receiver.
 29. Anapparatus as claimed in claim 16 and comprising a programmable dataprocessor.
 30. A storage medium containing a program for a dataprocessor of an apparatus as defined in claim
 29. 31. A method asclaimed in claim 9 and comprising the further step of processing theacquired seismic data using the de-signature and de-multiple operatorthereby to attenuate or remove seismic events arising from multiplereflections.
 32. A method as claimed in claim 9 and further comprisingthe steps of: emitting multi-component seismic energy at a sourcelocation; and acquiring seismic data at a multi-component seismicreceiver located at a greater depth than the source location.
 33. Amethod as claimed in claim 10 and further comprising the steps of:emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location.
 34. A method as claimed inclaim 12 and further comprising the steps of: emitting multi-componentseismic energy at a source location; and acquiring seismic data at amulti-component seismic receiver located at a greater depth than thesource location.
 35. A method as claimed in claim 13 and furthercomprising the steps of: emitting multi-component seismic energy at asource location; and acquiring seismic data at a multi-component seismicreceiver located at a greater depth than the source location.
 36. Amethod as claimed in claim 14 and further comprising the steps of:emitting multi-component seismic energy at a source location; andacquiring seismic data at a multi-component seismic receiver located ata greater depth than the source location.
 37. A method as claimed inclaim 10 and comprising the further step of processing the acquiredseismic data using the de-signature and de-multiple operator thereby toattenuate or remove seismic events arising from multiple reflections.